/* spectral_div_main.c */

#include "../common/spectral_div.h"

#ifdef MTRACE
#include <mcheck.h>
#endif

int main(int argc, char** argv)
{
	int result = 0, i = 0, g1_size = 0, g2_size = 0, g1_i = 0, g2_i = 0;
	int return_value = 1, elements_no = 0;
	sparse_matrix_lst* A = NULL, *A_g = NULL;
	int *g = NULL, *g1 = NULL, *g2 = NULL;
	elem *s = NULL, *k = NULL, *k_g = NULL, M = 0; 
	elem **B = NULL, **gen_B_g = NULL;
	elem *f = NULL, *res_vec = NULL;
	elem vec_mat_vec_res = 0;
	elem percision = 0;
	eigenpair EP;
	
	if (argc != 4)
	{
		result = printf("\n usage: spectral_div "
						"<adjacency-mat-file> <group-file> <precision>\n");
		assert(result >= 0);
		return -1;
	}
        
        #ifdef MTRACE
        mtrace(); /* Starts the recording of memory allocations and releases */
        #endif

	A = read_adjacency_matrix(argv[1]);
	if (A == NULL)
	{
		return_value = -1;
		goto main_cleanup;
	}

	g = (int*) calloc(A->n, sizeof(int));
	assert(g != NULL);

	elements_no = read_group_file(argv[2], g);

	if (elements_no == -1)
	{
		return_value = -1;
		goto main_cleanup;
	}

	B = create_modularity_mat(A);
	k = calc_nodes_degrees(A);
	k_g = create_k_g(k, g, elements_no);
	M = count_edges(k, A->n);
	f = create_f(B, g, A->n, elements_no);
	gen_B_g = create_generalized_modularity_matrix(B, g, elements_no, f, A->n);
	A_g = create_A_g(A, g, elements_no);

    percision = atof(argv[3]);
	EP = find_eigenpair(A_g, gen_B_g, k_g, f, M, g, elements_no, percision);

	if (IS_POSITIVE(EP.eigenvalue))
	{
		s = (elem*) calloc(elements_no, sizeof(elem));
		assert(s != NULL);

		compute_s(s, g, EP.eigenvector, elements_no);

		div_indices(s, g, &g1, &g2, &g1_size, &g2_size, elements_no);

		res_vec = (elem*)calloc(elements_no, sizeof(elem));
		assert(res_vec != NULL);

		mult_vector_matrix(res_vec, s, gen_B_g, elements_no);

		vec_mat_vec_res = mult_vector_vector(res_vec, s, elements_no)/2;
		free(res_vec);
	
		result = printf("%lf\n", vec_mat_vec_res);
		assert(result >= 0);

		print_groups(g1, g2, g1_size, g2_size);
	}
	else
	{
		result = printf("0\n");
		assert(result >= 0);

		print_group(g, elements_no);
	}

main_cleanup:
	if (g != NULL)
		free(g);
	if (B != NULL)
		free_matrix(B, A->n);
	if (f != NULL)
		free(f);
	if (gen_B_g != NULL)
		free_matrix(gen_B_g, elements_no);
	if (g1 != NULL)
		free(g1);
	if (g2 != NULL)
		free(g2);
	if (A != NULL)
		free_sparse_matrix_lst(A);
	if (k != NULL)
		free(k);
	if (k_g != NULL)
		free(k_g);
	if (EP.eigenvector != NULL)
		free(EP.eigenvector);

        #ifdef MTRACE
        /* End of memory allocation checks. */
        muntrace();
        #endif

	return return_value;
}
